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j  FREQUENCY  WAVENUMBER  BEAMFORMING  BY  USE  OF  THE 
TWO  DIMENSIONAL  FOURIER  TRANSFORM  , 

/ 

D.A.  Gray 


SUMMARY 

The  concept  of  frequency  wavenumber  beamforming  of  a  linear 
array  of  equi spaced  receivers  is  reviewed.  The  frequency 
wavenumber  spectrum  can  be  evaluated  by  using  a  two- 
dimensional  Fast  Fourier  Transform  and  the  properties  of 
the  resulting  estimate  of  the  frequency  wavenumber  spectrum 
are  described.  Particular  attention  is  given  to  the 
distortion  of  the  spectral  power  estimates  due  to  leakage 
of  negative  frequency  components. 
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1.  INTRODUCTION 

Beamforming  may  be  described  loosely  as  the  coherent  processing  of  an  array  of 
receivers  in  order  to  ascertain  the  arrival  directions  of  the  incident  signals 
or  to  increase  the  signal -to-noise  ratio  of  a  signal  from  a  given  direction. 

In  time  domain  beamforming  this  is  achieved  by  time  delaying  the  receiver  outputs, 
so  that  they  are  brought  into  phase  coherence  for  a  plane  wave  arriving  from  a 
given  direction,  and  subsequently  summing  them.  If  the  receiver  outputs  are 
first  narrowband  filtered  then  this  time  delay  may  be  effected  by  a  phase  multi¬ 
plication.  This  process  is  termed  frequency  domain  beamforming. 

However  the  choice  of  phase  factors  in  frequency  domain  beamforming  need  not  be 

restricted  to  those  corresponding  to  real  angles  of  arrival.  As  discussed  by 

Booker  and  Clemmow(ref. 1,2,3)  the  variable  wavevector  is  often  used  instead  of 

angular  coordinates.  Thus  the  total  power  incident  on  an  array  may  be  estimated 

as  a  function  of  frequency,  f,  and  the  wavevector  components  k  ,  k  and  k  (see 

x  y  z 

references  4  and  5) .  It  should  be  noted  that  only  three  of  the  above  four 
variables  are  independant.  This  power  estimate  is  teamed  the  frequency  wave- 
number  power  spectrum  (although  to  be  pedantic  it  should  be  teamed  the  frequency 
wavevector  power  spectrum) . 

For  a  linear  array  only  a  two  dimensional  foomulation  need  be  used.  Thus  the 

power  spectrum  is  estimated  as  a  function  of  the  variables  f,  k  and  k  .  Only 

x  y 

the  independent  variables  f  and  kx  (s  -p)  need  be  used.  In  this  paper  the 
liberty  of  referring  to  kx  as  wavenumber  rather  than  'the  component  of  the  wave- 
vector  along  the  axis  of  array*  is  taken  although  it  should  be  borne  in  mind  that 
wavenumber  is  strictly  defined*as 

From  figure  1  k  may  be  interpreted  as  a  spatial  frequency  and  for  a  linear  array 
of  equispaced  receivers  the  decomposition  into  an  angular  spectrum  of  plane  wave 
components  is  mathematically  equivalent  to  the  Fourier  analysis  of  a  sampled  time 
series  into  a  spectrum  of  single  frequency  components. 

All  phase  delays  correspond  to  particular  values  of  wavenumber  which  then  must 
be  interpreted  in  teams  of  the  physical  quantities  of  interest,  viz.,  angle, 
frequency,  and  velocity  of  propagation.  If  these  phase  delays  are  effected  at 
a  number  of  frequencies  then  the  resulting  power  spectrum  is  termed  the 
frequency  wavenumber  power  spectrum. 

For  a  linear  array  of  equispaced  receivers  it  is  shown  in  Section  2  how  the 
frequency  wavenumber  spectrum  can  be  written  as  a  two-dimensional  Fourier  trans¬ 
form  of  the  matrix  of  receiver  outputs.  For  appropriate  choices  of  the  number 
of  receivers  and  data  points  the  Fourier  transforms  may  be  effected  by  use  of  the 
Fast  Fourier  Transform  algorithm.  The  properties  of  the  frequency  wavenumber 
spectrum  are  discussed  in  Section  3  with  particular  attention  being  paid  to  their 
interpretation  in  terms  of  beamforming. 

In  conventional  time  series  analysis,  sidelobes  from  negative  frequency  components 
may  distort  low  frequency  components  in  the  estimate  of  the  power  spectrum  (see 
references  6  and  7) .  In  a  similar  way  the  negative  frequencies  introduced  by 
the  restriction  to  real  time  series  may  distort  other  low  frequency  components 
in  the  estimate  of  the  frequency  wavenumber  spectrum.  Some  measures  of  the 
sidelobe  level  and  the  power  distortion  due  to  this  leakage  are  derived  in 


Footnote:  Alternative  derivations  use  a  2ir  factor. 
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Section  4.  Finally,  in  Section  5  some  miscellaneous  remarks  on  the  treatment 
of  statistical  data,  the  derivation  of  the  spatial  autocorrelation  function  and 
methods  of  displaying  the  frequency  wavenumber  spectrum  are  made. 


2.  THE  FREQUENCY  WAVENUMBER  TRANSFORM 


2.1  Beamforming  fundamentals 


In  time  domain  beamforming  the  constructive  interference  of  signals  is 
effected  by  delaying  the  outputs  of  given  receivers  so  as  to  bring  into  phase 
coherence  any  signals  arriving  from  a  given  direction.  For  the  line  array 
of  N  equispaced  receivers  as  illustrated  in  figure  1  the  output  AC® ,  t)  of  a 
beam  'steered'  in  direction  0  is  given  by 


AC«,t) 


N-l 

I 


j=0 


x.  Ct  -  jAr) 


where  Ar  =  -s^ ,  d  is  the  separation  of  adjacent  receivers  and  c  is  the  prop¬ 
agation  velocity  of  the  signal  in  the  medium. 


In  order  to  determine  the  frequency  content  of  a  signal  incident  from  direction 
6  it  is  necessary  to  Fourier  analyse  A(0,  t) .  For  a  continuous  signal  Xt®,  f) 
an  estimate  of  the  complex  amplitude  of  a  signal  of  frequency  f  arriving  from 
a  direction  6 ,  is  given  by 

T 

A(S  ,t)exp(-2irift)  dt  (1) 


x(fl,f)  = 


iim 


T+oo 


2T 


lim 

T+oo 


(t  -  jAr)exp(-2irift) 


dt. 


The  power  at  this  frequency  and  direction  S(0,f)  is  then  given  by 

S(«,f)  =  AC®  ,  f)  A*  (0  ,  f) 

Similarly  the  total  power  incident  from  direction  0  is  given  by 


P(0) 


lim 

■p»oo 


A(0,t)  A* (0 , t) dt 
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2.2  Frequency  domain  beamforming 


It  is  well  known  that  for  a  narrow-band  signal  a  time  delay  can  be  effected 
as  a  phase  multiplication.  By  interchanging  the  summation  and  integral  in 
equation  (lj  this  principle  may  be  applied  to  the  time  delays  required  for 
beamforming.  It  then  follows  that 


N-l 


r  T 


X(d  : 

,f)  = 

lim 

'Y+oo 

N-l 

V 

) 

j=0 

N-l 

v 

L 

j=0 

the  wave 

-number 

,  is 

(lim  2^"  /  xj  (t)exp-(2Jrift)  dt)  exp (- 2m if  j At)  (2) 

T+oo  ^ -T 


Xj(f)  exp(-2irifjAr) 


k  -  f 


sind  sind 


(3) 


then 


N-l 


X(k,f)  *  \  Xj  (f)exp(-2>rijkd) 


(4) 


which  is  just  a  Fourier  series  of  the  Xj(f)'s.  The  new  variable,  the  wave- 

number,  has  a  direct  physical  significance  in  that  it  is  (apart  from  a  factor 
c)  the  spatial  frequency  of  the  incident  plane  wave  in  the  direction  of  the 
array  axis.  This  is  diagrammatically  illustrated  in  figure  1.  This 
interpretation  of  wavenumber  as  a  spatial  frequency  allows  a  very  useful 
analogy  between  time  series  analysis  and  frequency  domain  beamforming.  The 
spatial  series  X^ (f)( j  =  0,1... N-l  may  be  considered  as  equivalent  to  a  time 

series  of  N  data  points  and  the  resolution  into  spatial  frequencies  or  wave- 
number  expressed  in  equation  (4)  of  this  spatial  series  is  equivalent  to  the 
resolution  of  a  time  series  into  its  frequency  components.  This  idea  was 
first  put  forward  by  Booker  and  Clemmow(ref . 1) . 


FOOTNOTE:  As  discussed  in  the  introduction,  wavenumber  is  the  magnitude  of  a 
vector,  the  wavevector,  the  direction  of  which  is  orthogonal  to  the  wavefront. 
However,  since  only  one  dimension  is  considered  here  wavenumber  is  defined  in 
terms  of  the  projection  along  the  array  axis.  Alternative  definitions  involve 
a  2jt  factor. 
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2.3  Application  to  discrete  finite  time  series 

The  derivation  of  frequency  domain  beamforming  in  the  preceeding  section 
assumed  that  output  of  each  receiver  was  a  continuous  time  series  running 
from  t  =  -°°  to  t  =  +°°.  In  practice  the  time  series  is  only  observed  over  a 
finite  interval  and  often  only  sampled  data  is  available.  Under  these 
conditions  the  phase  delays  effect  not  a  true  but  a  cyclic  time  delay;  this 
is  demonstrated  in  Appendix  I.  As  discussed  in  Appendix  I  and  references  8 
and  9  the  error  associated  with  this  approximation  to  the  true  time  delays 
can  be  minimized  by  choosing  M  (the  number  of  data  points)  to  satisfy 


M 

c 


For  discrete  sampled  data  the  continuously  variable  vector  x^(t)  is  replaced 

by  the  matrix  X  =  x.(£t  )  abbreviated  to  x.;  where  j,  as  before,  runs  over 
jo  jt 

the  receivers  0,1,. 


,N-1  and  l  representing  the  time  samples  will  in  general 

The  correspond¬ 
ing  expression  from  X(k,f)  now  becomes 


take  the  values  0,1,..., M-l  and  is  the  sampling  interval. 


X(k,f) 


N-l  M-l 


j=0  1=0 


exp(-2ffikjd)  exp(-2mif£TQ)  . 


(5) 


While  X(k,f)  is  still  a  continuous  function  of  k  and  f  it  follows,  since  the 
transform  is  discrete  that  not  all  values  of  X(k,f)  are  independent.  The 
maximum  number  of  independent  values  can  be  generated  by  finding  a  complete 
set  of  orthogonal  vectors.  The  most  convenient  set  of  vectors  to  choose  are 
the  M  sine  and  cosine  waves  whose  periods  are  integer  subdivisions  of  the 
2 

observation  interval  and  the  N  complex  spatial  frequencies  whose  spatial 
periods  are  integer  subdivisions  of  the  array  aperture.  Once  these  orthogonal 
set  of  vectors  have  been  determined  the  value  of  X(k,f)  for  an  arbitrary  k  and 
f  can  be  expressed  as  a  linear  combination  of  these  orthogonal  ones.  The 
fundamentals,  Af  and  Ak  of  these  orthogonal  sets  of  frequencies  and  wavenumbers 
are  given  by 


and 


Akd 


1_ 

N  * 


Multiples  of  Af  and  Ak  arc  termed  the  frequency  and  wavenumber  bins  respectiv¬ 
ely.  The  set  of  MN  orthogonal  (in  a  two-dimensional  sense)  vectors  are 


for 


exp(-2irijij/N)  exp(-2ffi  vl^ 

H  =  1,2,. ..,N 


5 


WSRL-0162-TR 


and 


v  =  1,2, . . . ,M  . 


Denoting  X(Mk,  Mf)  as  X  ,  equation  (4)  reduces  to 

fAV 


N-l  M-l 


x.k  exp(-27riMj/N)  exp(-2ffiyd  /M) 


j  =0  k=0 


C7) 


This  can  now  be  readily  recognised  as  a  two-dimensional  discrete  Fourier 
transform.  The  power  spectral  density  S  is  given  by 

^  =  V  X£, 

and  S  v  Af  Ak  represents  the  power  incident  on  the  array  within  a  frequency 
range  of  (v  ±  h)  Af  and  a  wavenumber  range  of  (ju  ±  y  Ak. 


By  choosing  M  and  N  to  be  highly  composite  (eg  powers  of  two) ,  the  above 
transformation  may  be  effected  by  use  of  the  Fast  Fourier  Transform  algorithm 
and  so  the  number  of  operations  is  reduced  from  the  order  of  NJ M2  to  NM  logiM 
log2N. 

2.4  Physical  interpretation 

As  indicated  in  the  previous  section  S  is  an  estimate  of  the  power  density 
incident  on  the  array  corresponding  toM a  frequency  vAf  and  a  wavenumber  pAk . 
Since,  in  general,  the  physical  quantities  of  interest  are  frequency  and 
bearing,  then  wavenumber  should  be  reinterpreted  in  terms  of  these  parameters. 
From  equation  (3)  the  bearing  0 ,  corresponding  to  a  given  frequency,  f,  and 
wavenumber,  k,  is  given  by 


sin®  =  ck/f 


(8) 


Furthermore  denoting  B  as  the  angle  corresponding  to  the  v  -th  frequency  bin 
and  the  M  -th  wavenumber  bin  then 


sin  0 


\i.V 


cMk 

vAf 


cmNAt 


(9) 


There  are  two  interesting  consequences  that  follow  immediately  from  the  above 
equation. 
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(a)  Lines  of  constant  bearing 

From  equation  (8)  a  line  of  constant  bearing  is  represented  by  a 
linear  relationship  between  k  and  f  the  gradient  of  which  is  c/sin0  *. 
If  f  and  k  are  plotted  as  the  two  axes  of  a  two-dimensional  plot  then 
as  indicated  in  figure  2  lines  of  constant  bearing  will  run  radially 
outwards  from  the  origin  with  a  slope  proportional  to  c/sin0.  Since 
the  input  time  series  is  real  it  follows  that  the  negative  frequency 
components  are  mirror  images  of  the  positive  ones  with  opposite  k  (see 
later) .  In  figure  2  the  lines  corresponding  to  the  two  endfire  beams 
and  the  broadside  beam  have  been  shown.  The  frequency  f,  is  that 

'a 

corresponding  to  the  half -wavelength  of  the  array. 


(b)  Physical  region 

For  real  values  of  9,  sin  6  lies  between  ±1  but  as  can  be  seen  from 

ok 

equation  (3)  for  low  values  of  f,  l  .-j-l  can  be  greater  than 

unity,  eg,  f  =  0,  which  would  imply  that  I  sin  6  I  >  1.  In  the 
remainder  of  this  article  the  region  in  which  I  sin  6  |  <L,  ie, 

I  f  !  <  c.l  k  I  will  be  termed  the  "physical  region"  and  the  region 
where  I  f  I  >  c.  I  k  I  ,  the  "non-physical  region".  The  physical 
region  represents  plane  wave  signals  incident  on  the  array  from 
angles  lying  between  ±90  . 

The  crosshatched  areas  in  figure  2  represents  the  physical  region  and 
is  bounded  by  the  lines  corresponding  to  the  two  endfire  directions 
of  ±90°  .  The  number  of  frequency-wavenumber  bins  lying  within  the 
physical  region  can  easily  be  related  to  the  half -wavelength  of  the 
array.  At  the  half-wavelength  of  the  array 

\  '  2d  ’ 


and  substituting  in  equation  (9)  this  implies  that 


;in0 


N/2 


,P 


c£ 

f.  Md 

h 


_  M 

M 


for  P 


It  trivially  follows  that  since  I  p  I  <  j  all  the  wavenumber  bins 
correspond  to  physical  directions.  For  frequencies  f  less  than  f, 

Mf  * 

it  directly  follows  that  any  value  of  p  such  that  I  p  I  <  will 
correspond  to  a  real  direction. 


* FOOTNOTE :  In  the  figure  f  and  k  are  treated  as  continuous  variables;  for 
sampled  data  a  matrix  of  data  should  be  superimposed  on  the  diagram  as  shown. 
The  interdot  spacings  are  Ak  and  Af  along  the  wavenumber  and  frequency  axis 
respectively. 
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(c)  Non-physical  region 

As  discussed  above,  contributions  to  the  power  spectrum  lying  within 
the  "nonphysical  region"  do  not  correspond  to  a  real  signal  incident 
upon  the  array.  Some  possible  interpretations  of  components  in  this 
region  are: 

(i)  As  illustrated  in  figure  3  and  discussed  in  greater  detail  in 
Section  3.2,  leakage  from  either  the  frequency-time  transform 
or  the  wavenumber-space  transform  gives  rise  to  contributions 
within  this  region.  This  is  especially  true  of  the  negative 
frequency  sidelobes  which  are  discussed  in  greater  detail  in 
Section  4.  Sidelobes  of  the  frequency- time  transformation 
can  be  differentiated  from  lines  of  constant  bearing  since  they 
are  parallel  to  the  frequency  axis  and  do  not  run  radially  out¬ 
wards  from  k  =  0  as  do  lines  of  constant  bearing.  (See  later). 

(ii)  Often  effects  can  be  transmitted  through  the  array  structure  at 
velocities  differing  from  those  in  the  medium.  For  example, 
at  a  quarter  wavelength  the  physical  region  would  be  doubled  if 
the  velocity  of  propagation  through  the  array  was  half  that 
through  the  surrounding  medium.  This  is  particularly  common 
in  acoustic  arrays,  whereby  vibrations  of  the  array  structure 
can  propagate  along  the  array  and  excite  the  hydrophone 
receivers . 


3.  APPLICATION  OF  THE  TWO-DIMENSIONAL  FAST  FOURIER 
TRANSFORM  TO  BEAMFORMING 

In  this  section  some  properties  of  the  two-dimensional  discrete  transform  will 
be  discussed  with  particular  emphasis  being  placed  on  the  interpretation  of  these 
properties  from  a  beamforming  viewpoint . 

3.1  Symmetries 

Since 

N-l  M-l 

Sp„  =  I^T  ^Txjk  exp(-2*iPj/N)  •  exp(-27Tirk/M)  I1 
j=0  k=0 

it  immediately  follows  that 

Sp  ±  pN  ±  qM  =  Spp 

for  p  and  q  integers.  This  periodicity  in  frequency  and  wavenumber  is  a 
consequence  of  the  discrete  temporal  and  spatial  sampling  respectively.  This 
results  in  a  block  periodic  structure  as  indicated  in  figure  2.  Also  the 
x^j's  are  real  and  so  as  a  result  (see  Appendix  II)  it  follows  that 

SN-p  M-p  =  Spv 


SN-pt>  =  Sp  M-i'  • 


and 


(10) 
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The  physical  interpretation  of  equations (10)  is  that  the  positive  frequency 
components  from  a  direction  </>  and  negative  frequency  components  from  ^  +  180° 
are  identical  for  a  given  value  of  f.  At  a  given  frequency  the  positive  and 
the  corresponding  negative  wavenumber  terms  have  no  such  symmetry  since  they 
correspond  to  opposite  velocities  of  the  projected  wavenumber  along  the  array 
axis . 

Note  that  for  the  wavenumbers  to  run  contiguously  from  angles  corresponding 
to  to  |  the  terms  in  equation  (10)  should  be  rearranged  to 


and 


S '  N/2+pv  ~  V 


pv 


=  S. 


N/2+pv 


for  p  =  0,1, • • • ,N/2-l  . 


As  shown  in  Appendix  II  the  d.c.  and  Nyquist  terms  require  special  treatment. 

3.2  Response  to  a  complex  sinusoid 

As  discussed  in  Section  2.3  the  two-dimensional  Fourier  transform  is  a 
decomposition  of  the  matrix  X  into  orthogonal  plane  waves  which  are  harmonics 
of  exp(-27riAkx)  and  exp(-2friAft)  .  It  follows  that  any  plane  wave  input  which 
is  the  (k,p)th  harmonic  of  this  fundamental  will  contribute  zero  to  terms  in 
S  for  p  k  and  v  p.  For  plane  waves  which  are  not  exact  multiples  of 

r  ' 

the  spatial  and  temporal  harmonics  the  power  will  be  distributed  among  the 
Spv's  and,  by  analogy  with  time  series  analysis,  this  'leakage'  may  be  specified 

by  a  response  function.  This  response  function  Hp(>  (k,f)  is  the  S^resulting 

from  a  plane  wave  of  frequency  f  and  a  wavenumber  k.  The  response  function 
may  be  thought  of  as  a  generalization  of  the  polar  diagram  of  the  array  since 
it  includes  all  frequencies  and  extends  the  response  into  the  'non-physical' 
region  (see  Section  2.4).  It  also  holds  that  the  observed  output  from  an 
arbitrary  signal  is  given  by  the  convultion  of  the  true  frequency  wavenumber 
spectrum  of  the  signal  with  the  response  function. 

Furthermore,  since  for  a  complex  sinusoid  of  unit  amplitude  the  data  x.  is 
given  by  the  product 


exp(2irif jr  )  exp(2irikpd) 


it  follows  that  Hp(/ 


(f,k)  decomposes  into  the  product 


where 


Wp(>0  F^f) 


sin2  7rM(fr  -  jj) 

O  M 

M2  sin2  *(f tq  -  g) 


Wp(k) 


sin2irN(kd  -  £) 

N2  sin2  »r(kd  -  jj) 


(11) 


and 
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As  a  result  of  the  factorization  in  expression  (11)  many  of  the  properties  of 
the  response  function  may  be  treated  separately  in  spatial  and  temperal  co¬ 
ordinates. 


3.2.1  Leakage 

Leakage  can  occur  from  either  the  frequency  or  the  wavenumber  transform 
and  the  combined  effect  is  given  by  equation  (11). 


In  particular,  each  term  in  equation  (11)  is  of  the  form  G(^)  = 


sin 


P  (0-*o) 


Pas'in2  ($-^q) 


where  P  is  either  the  number  of  hydrophones  or  time  samples  and  <t>  represents 
either  2mk  or  2irf.  This  function  is  plotted  in  figure  4  where  A0,  the 
separation  of  the  zeros  and  also  the  bin  width,  is  either  the  reciprocal  of 
the  observation  interval  or  the  array  aperture.  The  leakage  into  the 
various  bins  is  determined  by  the  value  of  the  appropriate  response  function 
at  the  multiples  of  .  For  example  in  figure  4  when  0O(1)  is  an  exact 


multiple  of  A0  there  is  no  leakage  into  other  bins  since  they  correspond  to 
the  zeros  of  G(tf) .  The  worst  possible  case  is  when  $D( 2)  lies  halfway 


between  two  bins  as  illustrated  in  figure  4  and  the  power  is  equally 
distributed  between  the  two  adjacent  bins.  In  each  case  the  power  is 
approximately  3.9  dB  down  on  the  true  power.  In  the  extreme  case  of  both 
the  frequency  and  wavenumber  lying  exactly  between  adjacent  cells  it  follows 
from  expression  (11)  that  the  power  in  the  adjacent  cells  is  reduced  by 
7.8  dB. 


There  is  an  interesting  consequence  of  the  multiplicative  effect  of 
equation  (11)  which  has  a  simple  physical  interpretation.  The  polar 
diagram  as  a  function  of  k  is  only  dependent  on  the  frequency  f  through  a 
scaling  factor  which  is  just  F (f) .  The  physical  interpretation  of  this 
is  that  frequency  components  which  'leak'  from  the  time -frequency  trans¬ 
formation  are  beamformed  in  k-space  but  with  a  reduced  amplitude.  The 
invariance  of  the  response  function  W(k)  on  frequency  is  one  of  the 
advantages  of  using  a  frequency  wavenumber  as  opposed  to  a  frequency  bear¬ 
ing  representation. 

To  illustrate  these  considerations,  two  response  functions  for  a  32  element 
array  are  plotted  in  figures  3  and  5.  These  functions  are  plotted  as  a 
continuous  function  of  wavenumber  but  the  frequencies  were  evaluated  at 
multiples  of  the  frequency  bin  Af  where  the  highest  frequency  corresponded 
to  X/2  of  the  array.  In  figure  3  the  complex  sinusoid  was  an  integer 
multiple  of  Af  and  no  frequency  leakage  is  apparent.  However  in  figure  5 
since  the  frequency  was  chosen  to  lie  halfway  between  two  adjacent  bins, 
leakage  caused  significant  contributions  to  the  resulting  spectrum  at  all 
frequencies.  The  radical  difference  between  these  two  spectra  show  the 
importance  of  choosing  a  frequency  resolution  which  is  finer  than  that  of 
the  expected  line  width  of  data  or  weighting  the  time  series  to  reduce 
leakage. 

Another  advantage  of  using  the  frequency  wavenumber  representation  is 
shown  in  figure  5.  If  the  frequency  resolution  is  inadequate  and  leakage 
does  occur  then  peaks  of  the  ' leaked' frequency  terms  run  parallel  to  the 
frequency  axis.  This,  even  with  noise  present,  serves  to  distinguish 
them  from  broadband  targets  of  a  constant  bearing  which  would  run  radially 
from  the  centre,  (see  also  figure  6). 
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3.2.2  Ripple 

Often  it  is  required  to  estimate  more  than  either  the  M/2  frequency  bins 
or  the  N  wavenumber  bins.  One  important  reason  for  this  is  that  due  to 
the  fact  that  the  spectra  corresponding  to  sinusoids  with  wavenumber  or 
frequency  components  lying  between  bins  will,  due  to  leakage,  have  a  ripple, 
the  depth  of  which  can  be  up  to  7.8  dB.  This  ripple  may  lead  to  false 
conclusions.  Another  reason,  particularly  valid  for  the  spatial  transform, 
is  that  the  number  of  spectrum  points  is  insufficient  and  it  may  be 
desirable  to  interpolate  between  them.  A  simple  solution  to  either 
problem  is  to  double  the  number  of  beams  steered;  this  reduces  the  spatial 
ripple  to  less  than  1  dB. 

Fortunately  a  simple  way  of  generating  these  redundant  beams  is  to  add  a 
block  of  N  (or  M  in  the  time  domain)  zeros  to  the  appropriate  row  or 
column  of  the  matrix.  The  required  number  of  beams  are  generated  by  an 
FFT  or  length  2N  (or  2M) .  If  further  interpolation  is  required  then  more 

blocks  of  zeros  should  be  similarly  added. 

3.2.3  Spatial  aliasing 

It  is  often  the  case  that  the  sampling  rate  is  sufficiently  high  that  the 
frequency  transformation  will  produce  frequencies  corresponding  to  wave¬ 
lengths  shorter  than  the  half-wavelength  of  the  array.  It  is  well  known 
that  the  polar  diagram  of  an  array  at  these  higher  frequencies  may  contain 
extra  main  lobes.  By  analogy  with  time  series  analysis,  this  phenomenon 
is  termed  spatial  aliasing  and  once  again  has  a  convenient  interpretation 
in  the  wavenumber  representation.  As  discussed  in  Section  3.1  the 
frequency-wavenumber  spectrum  is  periodic  in  wavenumber  and  in  a  similar 
way  to  time  series  analysis  aliased  wavenumbers  fold  back  into  the  wave- 
number  spectrum.  However,  since  the  positive  and  negative  wavenumbers 
are  not  symmetric  as  in  time  series  analysis  the  aliased  peaks  of  positive 
wavenumber  will  be  present  as  negative  wavenumber  components  and  vice 
versa.  An  example  of  this  is  shown  in  figure  6  where  the  frequency  range 
of  a  32  element  array  has  been  extended  to  3^/4  introducing  aliased  positive 
wavenumbers  for  a  broadband  target. 

Once  again  the  frequency-wavenumber  representation  is  an  appropriate  mode 
for  this  analysis.  This  is  because  lines  of  constant  bearing  are  folded 
back  as  straight  lines  which  no  longer  run  radially  out  from  the  origin 
but  are  inclined  at  a  characreristic  angle.  This  is  illustrated  in 
figure  6  for  a  simulated  broadband  source. 

3.2.4  Windowing  and  shading 

As  can  be  seen  from  the  examples  of  Section  3.2.2  the  sidelobes  of  either 
the  frequency  or  wavenumber  transformation  may  often  be  too  high.  In 
particular  the  leakage  of  power  from  a  strong  signal  through  these  side- 
lobes  may  distort  and  even  hide  a  weaker  signal.  One  technique  to  reduce 
these  sidelobes  is  to  multiply  each  data  sample  by  a  weight  which  is 
designed  to  reduce  the  sidelobes.  This  is  termed  windowing  in  time  series 

analysis  and  shading  in  array  beamforming.  These  are  a  number  of  well 
known  weights  among  which  are  the  Hamming,  Parzen  and  Tchebyshev  which  can 
be  applied  in  either  the  temporal  or  spatial  domain.  All  of  these  although 
reducing  sidelobes,  have  the  inevitable  effect  of  increasing  the  width  of 
the  main  peak.  A  more  complete  description  of  the  effect  of  shading  is 
given  in  reference  11. 
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A  general  window  and  shading  scheme  can  be  represented  by  where 

and  to,  are  the  appropriate  temporal  and  spatial  weights  respectively. 

The  matrix  x^v  is  then  replaced  by  x . 

3.3  The  relation  of  the  wavenumber  bins  to  beamwidths 

As  discussed  in  Section  3.1  there  are  N  independent  wavenumber  beams  and  as 
shown  in  the  diagrammatic  representation  of  H(k)  in  figure  4  adjacent  beams 
overlap  each  other  at  their  approximate  3.9  dB  points.  These  beams  are 
equispaced  in  wavenumber  for  all  frequencies.  Since  the  size  of  the 
'physical  region'  varies  linearly  with  frequency  the  invariance  of  H(k)  has 
a  number  of  important  consequences  when  interpreted  in  terms  of  angles  and 
are  discussed  below. 

3.3.1  Frequency  broadening 

Frequency,  wavenumber  and  angle  are  related  by 


k  =  2m  iJinl 

c 

An  approximation,  particularly  valid  for  small  0  is  that  the  change  in 
angle  60  resulting  from  a  change  in  wavenumber  5k  is  given  by 

50  =  5kc/2Jrf  .  (12) 

If  6k  is  taken  to  be  the  constant  wavenumber  bin  then  it  follows  that  the 
corresponding  angular  separation  50,  or  beamwidth  as  it  is  commonly  known, 
will  increase  as  the  frequency  is  decreased  according  to  equation  (12). 

This  can  be  shown  to  be  true  for  all  0  and  is  an  expression  of  the  well- 
known  fact  that  the  angular  resolution  of  an  array  decreases  with  decreas¬ 
ing  frequency. 

3.3.2  Angular  broadening  at  a  given  frequency 

At  broadside  sin  0^0  and  so  the  beams  are  approximately  equispaced  in 
angle.  However  towards  endfire  an  increment  of  one  wavenumber  bin  will 
correspond  to  a  large  angular  variation.  Thus  the  beamwidth  of  the  end- 
fire  beams  will  be  broader  than  those  at  broadside.  Spacing  beams  such 
that  the  adjacent  beams  overlap  at  their  3.9  dB  points  is  a  simple  way  of 
ensuring  a  complete  coverage  with  a  minimum  of  redundancy  and  the  frequency- 
wavenumber  transform  provides  a  natural  formalism  for  attaining  this. 


4.  NEGATIVE  FREQUENCY  DISTORTION 

In  assessing  the  properties  of  the  frequency  wavenumber  spectrum,  frequent  use 
was  made  of  the  response  function,  defined  to  be  the  spectrum  corresponding  to 
a  complex  sinusoidal  input.  However  a  more  realistic  representation  of  the  data 
is  as  a  real  sinusoid. 

The  results  of  time  series  analysis  shows  that  in  general  the  effect  of  this 
restriction  is  slight  but  can  be  significant  in  the  treatment  of  low  frequencies. 
In  particular  since  the  real  sinusoid  is  the  sum  of  two  equal -amplitude  complex 
exponentials  of  positive  and  negative  frequencies,  leakage  of  the  negative  power 
peaks  will  distort  the  estimate  of  the  positive  frequencies.  Since  the  energy 
'leaked'  through  a  particular  sidelobe  decreases  with  its  distance  from  the  main 
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lobe,  low  frequency  components  are  more  likely  to  be  distorted  than  higher 
frequency  ones.  This  is  illustrated  in  figure  7.  Whether  the  interference  is 
constructive  or  destructive  will  depend  on  the  phase  of  the  sinusoid.  The 
effect  of  this  in  conventional  time  series  analysis  has  been  extensively  discussed 
in  references  1  and  2.  In  the  remainder  of  this  section  the  effect  of  these 
negative  frequencies  on  the  beamforming  will  be  discussed. 

4.1  Response  to  a  real  sinusoid 

Ignoring  phase  terms  the  real  data  matrix  due  to  a  cosinusoidal  input  can  be 
written 


x. 

H 


+ 


) 


where 

xt  =  exp(2ffi(f ' jrQ  -  k'  Id)) 

J  V 

and 


+  * 


For  this  section  it  is  more  convenient  to  work  with  the  continuous  tranform 
Y(k,f)  defined  by  equation  (S)  rather  than  just  the  independent  terms  of 
equation  (7) . 


Substituting  in  equation  (5)  the  above  expression  for  x 


jl 


it  follows  that 


X(k,f)  =  X  (k,f)  ♦  X'(k,f) 


where 


X  (k,f)  = 


M-l  N-l 

^  exp(2n-i(f'-f)jro)  exp(-2»ri(k'-k)^d). 

j=0  1=0 


and 


M-l  N-l 

X’(k.f)  =  pj^j  ^  ^exp(2ffi(f'+f) jTQ)  exp(-2iri(k'+k)*d). 

j=0  1=6 


Defining 


S  (k,f)  =  X  (k,f)  X  (k,f) 


it  can  readily  be  shown  that 
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rsirtTM(f-f ')r  ^ 2 
S  (k,f)  =  1 - & 


siiwrN(k-k  'Jd'j  2 


MsiiOT  (f-f  ')tq  j  (^sinJr(k-k')d 


-sir 

^Isi 


Similarly  S  (k,f)  defined  by  X  (k,f)  X  (k,f)  is  given  by 


S"(k,f) 


sinJTM(f+f '  )t  ^ 
Msinff  (f+f '  )r 


sin7rN(k+k')d-N  2 


Nsinff  (k+k’)d 


The  corresponding  response  to  x.  is  given  by 


S(k,f)  =  (X  (k,f)  +X'(k,f))  (X  Ck,f)  +X'  (k,f)) 
which  ignoring  the  cross  terms  reduces  to 

S(k,f)  a  S+(k,f)  +  S"(k,f). 


(13) 


Thus  the  amplitude  of  the  negative  frequency  leakage  into  the  positive 
frequency  region  is  approximately  determined  by  equation  (13)  evaluated  for 
f  >  0.  The  peak  of  S(k,f)  as  for  f  >  0  as  a  function  of  k  occurs  at  -k. 

This  affords  a  simple  physical  interpretation  that  the  negative  frequency 
components  which  leak  into  the  positive  frequency  region  are  beamformed  as 
normal  sinusoids  with  a  wavenumber  of  -k  i.e.  incident  upon  the  array  from 
the  reflected  angle  about  broadside. 

To  qualitatively  illustrate  these  features  both  S+(k,f)  and  S(k,f)  are  plotted 
in  figures  5  and  8  respectively  for  an  array  of  32  receivers  and  16  independent 
frequencies  varying  from  0  to  X/2  of  the  array.  Both  of  these  functions  have 
only  been  evaluated  for  the  discrete  values  of  f  corresponding  to  centre  cell 
when  using  a  32  point  frequency  transform. 

4.2  Negative  frequency  artefacts 

The  height  of  these  negative  frequency  sidelobes  which  leak  into  the  positive 
frequency  region  (termed  artefacts)  will  depend  not  only  on  N  and  M  (the 
number  of  data  points  and  receivers)  but  also  on  the  frequency  f'  and  wave- 
number  k'  of  the  incident  sine  wave.  The  envelope  of  the  negative  frequency 
beams  is  given  by 


1  1  1 
4 (MN) 2  sinJff(f  +  f’)rQ  sinJff(k  +  k’)d . 


(14) 


A  useful  measure,  C,  of  this  leakage  is  given  by  the  ratio  of  the  maximum 
value  of  the  envelope  for  f  >  0  to  the  peak  of  S+(k,f)  ie 


C 


10  log/ (S  '  (k , 

U*(k, 


max 

max 
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In  most  cases  the  approximation  siitfrff  +  f')ro  —  *(.£  +  f') rQ  holds  and  hence 

the  maximum  of  expression  (14)  occurs  for  f  =  0  and  k  =  -k'.  Thus  C  can  be 
shown  to  reduce  to 


C  =  -20  log  Mfff'r  . 

o 


However 


Mr  =  T 
o 

where  T  is  the  observation  interval  and  so  f'T  is  just  the  number  of  periods 
within  the  observation  interval.  Thus 

C  =  -20  log  irNT 


where  is  the  number  of  periods  of  the  sine  wave  within  the  observation 
interval . 

The  extension  of  the  above  argument  to  the  case  of  the  discrete  transform 
where  only  independent  values  of  the  two-dimensional  transform  are  evaluated 
can  easily  be  effected  by  means  of  the  sampling  theorem.  For  an  input 
sinusoid  whose  frequency  is  centrecell  then  no  artefacts  will  occur  since  all 
other  frequencies  occur  at  zeros  of  the  response  function.  The  worst  case 
occurs  when  the  frequency  of  the  sinusoid  lies  exactly  halfway  between  bins 
and  expression  C  (with  an  extra  3.9  dB  subtracted  to  allow  for  the  equiparti- 
tioning  of  S+(k,f)  into  adjacent  cells)  describes  the  heights  of  the  negative 
frequency  artefacts. 

4.3  Power  distortion 


When  the  'leaked'  negative  frequencies  artefacts  are  beamformed  the  resulting 
sidelobes  will  interfere  with  the  main  positive  frequency  lobe  and  as  a  result 
will  distort  an  estimate  of  the  power.  A  measure,  P,  of  the  distortion  in 
the  peak  power  estimate  of  a  real  sinusoid  can  be  defined  to  be  the  ratio  of 
S‘(k,f)  to  S+(k,f)  when  both  evaluated  at  k  =  k'  and  f  =  f ' . 

Thus 


p  ■  10 


sin  2mMf'T 

0 

sin  2»rNk'd 

M  sin  27rf,T 

0 

N  sin  2>rk'd 

This  varies  from  zero  when  both  f'  and  k'  are  centre  cell  to  a  maximum  when 
both  are  exactly  between  cells.  In  the  latter  case,  P  is  given  by 


P  =  -20  log  M  sin  2nf'r  N  sin  2>rk'd 

6  o 


which  can  be  approximated  by 
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P  =  -20  log  4JT2  M^Nt 

where  NT  is  the  number  of  periods  of  the  sine  wave  in  the  observation  interval 
and  M^  is  the  number  of  wavelengths  projected  onto  the  array  aperture.  (This 
letter  quantity  obviously  varies  with  incident  angle) .  Thus 

P  =  °X  Ct  • 


From  this  equation  it  can  be  seen  that  the  distortion  will  be  greatest  at  low 
frequencies  and  angles  around  broadside. 


5.  MISCELLANEOUS  COMMENTS 


S.l  Statistical  data 


In  the  preceding  chapters  the  properties  of  the  frequency  wavenumber  spectrum 
have  been  discussed  with  respect  to  deterministic  signals.  In  order  to  treat 
random  data,  some  sort  of  averaging  is  necessary  to  reduce  the  variance  of 
the  power  estimates  which  are  now  random  variables.  Two  methods  for  effect¬ 
ing  this  are  possible. 


(a)  In  time  series  analysis  the  ratio  of  the  variance  to  the  mean  of 
the  power  spectral  estimate  is  reduced  by  segmenting  the  time 
series  into  a  number  of  blocks,  estimating  the  power  spectrum  for 
each  block  and  then  averaging  the  power  spectra.  In  a  similar 
way  the  time  series  from  each  receiver  are  segmented  into  blocks 
of  M  data  points  and  the  frequency  wavenumber  spectrum  is  estimated 
for  each  block.  The  resulting  blocks  are  averaged  and  the  averaged 
power  spectral  estimates  will  have  X2  distribution,  the  ratio  of  the 
variance  to  the  mean  of  these  power  estimates  will  be  reduced  by  a  factor 
of  >/~P  where  P  is  the  number  of  averages.  Since  the  array  aperture  is 
limited  in  space  the  segmenting  into  blocks  is  effected  in  the  time 
domain. 


An  attractive  feature  of  frequency  wavenumber  analysis  is  that  if 
only  independent  frequency  and  wavenumber  bins  are  used  then  detection 
statistics  may  readily  be  evaluated. 


(b)  Multiplying  the  data  matrix  by  a  window  w^j  can  easily  be  shown 
to  be  equivalent  to  convolving  the  Fourier  spectrum  of  the  X^'s  with 
the  Fourier  spectrum  of  the  W^'s.  Since  the  two-dimensional  con¬ 
volution  is  essentially  a  weighted  averaging  of  the  S^'s  over 

adjacent  bins  it  follows  that  the  variance  of  this  convolved  spectral 
estimate  will  be  reduced.  However,  the  bias  (or  lack  of  resolution) 
will  be  correspondingly  increased. 


5.2  Spatial  autocorrelation  function 

In  Section  2  it  was  shown  that  the  wavenumber  amplitude  spectrum  X(k,f)  is 
the  Fourier  transform  of  the  narrowband  aperture  distribution  of  signal 
Xj(f).  It  follows  directly  that  S(k,f)  *  X(k,f)  X*(k,f)  is  the  Fourier 

transform  of  the  function 


J 
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*  £ 
i=n 


x*Cf)  x.  ff) 

3  J+P 


(15) 


This  is  easily  seen  to  be  the  discrete  autocorrelation  function  r^(f)  of  the 


X.(f)'s. 


It  follows  that  the  frequency  wavenumber  spectrum  could  be 


alternatively  derived  as  the  Fourier  transform  of  the  narrowband  spatial  auto¬ 
correlation  function  at  each  frequency  of  interest.  If  an  assumption  of 
spatial  stationarity  is  made,  an  alternative  definition  of  r  (f)  would  be 


rp(f)  =  <X}(f)  X.+p(f)  >  (16) 

where  O  denotes  the  ensemble  average. 

Apart  from  the  Bartlett  window  in  equation  (15)  in  the  ensemble  limit  the 
estimators  given  by  (15)  and  (16)  are  equivalent.  In  practice  however  the 
ensemble  average  will  usually  be  replaced  by  a  time  average  over  a  finite 
time  interval  and  so  expressions  (15)  and  (16)  will  no  longer  be  equivalent. 
This  has  an  effect  on  a  number  of  non-linear  techniques  for  estimating  the 
frequency  wavenumber  spectrum  which  are  based  on  the  use  of  the  spatial  auto¬ 
correlation  function  (15)  or  the  cross  power  spectral  matrix. 

5.3  Displays 

5.3.1  Choice  of  axes 

A  common  alternative  to  frequency  wavenumber  analysis  is  to  only  form 
beams  corresponding  to  the  physical  region  and  to  form  these  beams  at 
equal  angular  or  sin#  intervals.  A  lengthy  comparison  of  the  different 
displays  is  given  in  reference  10  where  the  frequency  wavenumber  plots  are 
referred  to  as  a  distortionless  Bearing  Frequency  display.  An  advantage 
of  the  frequency  wavenumber  display  is  that  a  broadband  target  will  appear 
as  a  straight  line  whereas  on  a  frequency  bearing  display  it  appears  as  a 
funnel  shape.  For  detection  purposes  the  eyeball  matched  filtering 
required  is  considerably  simplified. 

To  convert  a  frequency  wavenumber  spectrum  into  a  frequency  -  sin0  spectrum, 
a  simple  interpolation  procedure  may  be  used:  For  each  frequency  a  number 

of  wavenunber  points  are  Fourier  Transformed,  zeros  are  augmented  to  this 
transform  and  then  the  inverse  transform  is  effected.  The  number  of  wave- 
number  points  used  in  the  forward  transform  should  at  least  cover  the 
physical  region  and  due  to  sidelobe  leakage  there  is  probably  considerable 
advantage  in  taking  more.  The  number  of  zeros  augmented  will  depend  on 
the  degree  of  interpolation  required.  P  computational  advantage  of  this 
interpolation  is  that  it  should  only  need  to  be  effected  at  the  end  of  an 
integration  period  on  the  averaged  spectrum. 

5.3.2  Display  technique 

The  method  chosen  to  display  the  intensity  as  a  function  of  the  variables 
frequency  and  wavenumber  has  been  to  use  a  'hide'  routine.  The  perspective 
of  this  can  be  varied  in  order  to  display  features  which  may  have  been 
blanked  out  by  strong  targets.  Two  alternative  methods  of  displaying  a 
function  of  two  variables  are  the  two-dimensional  intensity  modulated 
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displays  (or  LOFARgrams)  or  contour  plots  ('fraz')  or  equi-level  contours. 

5.4  Analysis  of  real  data 

To  illustrate  some  of  the  features  of  frequency  wavenumber  beamforming  128 
data  points  from  each  of  32  hydrophones  of  a  linear  array  were  processed 
using  a  two  dimensional  Fast  Fourier  Transform.  Zeros  were  appended  to  the 
end  of  each  spatial  transform  resulting  in  128  wavenumber  bins.  The 
frequency  wavenumber  power  spectrum  when  averaged  over  300  integrations  is 
shown  in  figure  9  where  the  40  independent  frequencies  range  from  zero  to  that 
corresponding  to  the  half-wavelength  of  the  array.  Features  such  as  the 
radial  characteristics  of  broadband  sources,  spatial  aliasing  of  low  frequency 
array  self  noise  and  leakage  into  the  nonphysical  region  are  readily  discernible. 


6.  CONCLUSIONS 

Provided  care  is  taken  to  avoid  the  errors  in  approximating  delays  as  discussed 
Appendix  I  and  the  leakage  through  negative  frequency  sidelobes,  the  two- 
dimensional  Fast  Fourier  Transform  has  been  shown  to  be  a  natural  and  efficient 
way  of  estimating  the  two-dimensional  frequency  wavenumber  spectrum.  The 
spacing  of  beams  which  naturally  overlap  at  their  3.9  dB  points  and  the  fact  that 
broadband  targets  appear  as  straight  lines,  are  two  particular  advantages  of 
using  a  frequency  wavenumber  analysis.  For  certain  problems  such  as  array  self¬ 
noise  the  extension  of  the  wavenumber  representation  into  the  non-physical  region 
is  very  useful. 
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APPENDIX  I 

PHASE  DELAYS  AS  A  CYCLIC  TIME  SHIFT 


Let  X.(f) 
I 


denote  the  frequency  output  of  the 


.th 

1 


receiver; 


it  is  defined  by 


Xj(f) 


M'l 


x.  (lro)  exp(-27riflro) 


(1.1) 


where  is  the  sampling  interval. 

Suppose  it  is  required  to  steer  a  beam  in  a  direction  <t>  defined  by 

_  dsin <t> 


where  p  is  some  integer.  (The  restriction  to  a  steering  angle  such  that  the 
interelement  time  delay  is  an  integral  multiple  of  the  sampling  rate  is  not 
essential  but  simplifies  the  following  argument) . 

When  steering  in  this  direction  X.(f)  is  multiplied  by  a  phase  factor  of  the 
form  1 


exp(-2JTifpro) 

However  since  the  transform  in  equation  (1.1)  is  finite  the  phase  factor  corresponds 
not  to  a  true  time  delay  but  rather  a  cyclic  time  delay.  Thus  if  the  input  series 
used  in  (1.1)  is 


V  V  X2  .  *M-1 


where 


x.  denotes  x.fiT  ) 
l  j  o 


then  it  follows  that 


X.Cf)  exp(-2jrifpro) 


corresponds  to  a  cycled  version  i.e. 


th  .  . 
p  position 

1 


^-p  ^-p+l  ^-1  Xo  X1  .  ’Si-p-l 
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In  a  time  delay  beamformer  the  corresponding  series  would  be 

-p  -1  o  1  T4-p-l 

Thus  since  the  frequency  domain  beamformer  does  not  have  available  the  values 
x  p,  , . , ,  x  j  it  approximates  them  by  x^  ,  ....  x^ 

The  magnitude  of  this  error  depends  on  a  number  of  factors;  the  frequency  of 
the  signal,  the  arrival  angle  and  the  length  (MrQ)  of  the  data  record  used  to 

effect  the  transform  (1.1).  The  problem  has  been  addressed  in  references  3  and 
4  and  to  minimize  the  effect  of  the  cyclic  time  delay  it  is  reconmended  that  the 
data  accumulation  period  (ie  Mtq)  should  be  much  greater  than  the  time  for  a 

signal  to  transverse  the  array.  This  time  is  longest  for  beams  steered  around 
endfire  and  so  the  condition  can  be  expressed  as 

Mr  »  N  -  • 

o  c 
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APPENDIX  II 

SYMMETRY  PROPERTIES  OF  THE  TWO-DIMENSIONAL  SPECTRUM 

The  number  of  independent  terms  in  the  power  matrix  S  are  derived  and  the  correct 
treatment  of  the  d.c.  and  Nyquist  terms  is  given. 

Symmetries 

Defining  the  following  cosine  and  sine  transforms 


v  2npj  2nkv 

>  =  \  xjk  C0S  N  C0S  ~M> 

3  >  k 


.  v  .  2npj  .  2ffkJ' 

V  =  \  Xjk  Sin  N  Sin  — ’ 
1  »k 


v  2ffpi  .  2nkv 

V  =  Xjk  C0S  N  sin  W  ’ 
3  >k  J 


(II. 1) 


and 


it  follows  that 


,  v  2mp  j  2nkv 

APv  =  \  xjk sin  -r cos  -r* 

J  »k 


X  =  a  -b  +  i  (c  +  d  1. 
pv  pv  pv  pv  pvJ 


From  equations  (II. 1)  the  symmetries  such  as 


apv  apM-v 


V, 


■pv 


^-pM- 


can  easily  be  derived,  with  corresponding  relationships  holding  for  the  b  ,  c 

and  d  .  An  immediate  consequence  of  the  reality  of  the  x..  's  is  the  relation 
pv  3 K 


XN-pM-*' 


and  consequently  that 


SN-p  M-v  =  Spv 


It  also  holds  that  since  S  is  given  by 


S  =  (a  -  b  )■*  +  (c  ♦  d  )3 
pv  K  pv  pvJ  pv  pvJ 


(II. 2) 


N-pv 


=  S 


PM-v 


that 
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These  symmetries  combine  to  reduce  the  number  of  independent  elements  by  a  factor 
of  2. 


D.C.  and  Nyquist  terms 

The  d.c.  and  Nyquist  terms  (corresponding  to  subscripts  0,  ^  or  j)  are  usually 

calculated  as  the  real  and  imaginary  parts  of  the  Xqo  element.  To  estimate  the 

power  in  these  components  by  equation  (II. 2)  is  incorrect.  Assuming  that  the 
d.c.  and  Nyquist  components  are  stored  in  the  real  and  imaginary  parts  of  the 
zero-th  term  the  correct  expressions  for  the  powers  are 


and 


ov 

3N/2^ 

S 

VO 


a2  +  c2 
ov  ov 


b2  +  d2 
c*>  ov 


a2  ♦  d2 

VQ  1>Q 


v  t  0 


SfM/2 


J 


For  v  =  0  the  expressions  are  further  contracted  to 


and 


S  =  a2 
oo  oo 

SN/2  M/2  =  b°° 


As  a  result  of 
power  terms  in 


this 
MN 
2  + 


decomposition  and  the  symmetries,  the  number  of  independent 

1. 


i  - 


POWER  (dB) 
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Figure  3 


Figure  3.  Side  lobe  leakage  into  non  physical  region 
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Figure  S.  Response  function  of  a  maximally  frequency 
leaked  sinusoid 


Figure  6.  Spatial  aliasing 


i 


Figure  9.  Frequency  wavenumber  spectrum  for  real  data 
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